Energy transfer process in gas models of Lennard-Jones interactions 
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We perform simulations to investigate how the energy carried by a molecule transfers to others in 
an equilibrium gas model. For this purpose we consider a microcanonical ensemble of equilibrium 
gas systems, each of them contains a tagged molecule located at the same position initially. The 
ensuing transfer process of the energy initially carried by the tagged molecule is then exposed in 
terms of the ensemble-averaged energy density distribution. In both a 2D and a 3D gas model with 
Lennard-Jones interactions at room temperature, it is found that the energy carried by a molecule 
propagates in the gas ballistically, in clear contrast with the Gaussian diffusion widely assumed in 
previous studies. A possible scheme of experimental study of this issue is also proposed. 

PACS numbers: 05.60.Cd, 51.10.+y, 51.20.+d 
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One important task of statistical mechanics is to un- 
derstand various transport processes. A well-known suc- 
cessful example is the self-diffusion of gas molecules: Due 
to Einstein's 1905 work, it has been widely accepted 
that a particle in a gas undergoes the Brownian motion 
[H, which can be modeled essentially with the random 
walk 2] and the resulting probability distribution func- 
tion (PDF) follows the diffusion equation [l|, Q ■ Because 
of its fundamental importance for various scientific dis- 
ciplines, the study of Einstein's theory has never ceased. 
Very recently, the direct experimental measurement of 
the instantaneous velocity of a Brownian particle in a 
gas has been realized, and the random walk picture was 
confirmed with high precision 

Another question of fundamental importance is how 
the energy carried by a molecule transfers with time, 
which is a key step towards understanding the macro- 
scopic energy transport. In general, the existing theories, 
taking for example the Helfand theory approached 
this issue by simply extending the random walk picture 
of a Brownian particle, and predicted a Gaussian en- 
ergy density distribution as well. However, it should be 
pointed out that whether random walk is the underlying 
mechanism of the energy dispersion has never been ex- 
amined experimentally nor numerically in a gas or more 
generally in fluids. 

Unlike in the study of the self-diffusion (or mass dif- 
fusion) where the position of a particle can be traced 
accurately [IHl] (and now even its instantaneous veloc- 
ity can be measured Q), a key difficulty in the study of 
the energy transfer is that it is hard to trace the energy 
transferred from particle to particle. The mass diffusion 
can be explored by focusing on the trace of an individual 
particle, but by nature, the energy transfer is a collec- 
tive behavior involving all the molecules related by the 
transferred energy. 

In this work we perform an equilibrium molecular 
dynamics investigation to explore how the energy of a 
molecule may transfer in a gas. We will restrict our- 



selves to a 2D gas model, but it has been verified that 
in its 3D counterpart the results remain qualitatively the 
same. We assume that the gas is composed of only one 
kind of molecule with diameter a and mass m. The setup 
consists of a square of area S with periodic boundary con- 
ditions and N molecules moving inside. The interaction 
between molecules is given by the Lennard-Jones poten- 
tial and the Hamiltonian of the system reads 
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where e is a constant governing the interaction strength 
and rij denotes the distance between molecules i and j. 
In our calculations the dimensionless parameters are set 
tobe<7 = l,m = l,e = l and the Boltzmann con- 
stant fee = 1; In addition, the number of the molecules 
N = 2500 and the area S = 200 x 200 are adopted. An- 
other important parameter is the temperature, which is 
fixed at T = 2.5, a value that corresponds to the room 
temperature with other adopted parameter values. To 
make the simulations more efficient, the potential energy 
between two molecules is approximated by zero when 
their distance is larger than r c — 3.5, as convention- 
ally adopted in the molecular dynamics studies of gases. 
Given these, the evolution of a system can be simulated 
straightforwardly. In our calculations the 7-order Runge- 
Kutta algorithm [9] with step 0.01 is employed. 

The ensemble is prepared through the following three 
steps: (i) First, a "seed" equilibrium system of N 
molecules at temperature T is prepared by evolving a 
system for a long enough time (> 1 x 10 6 ) from a prop- 
erly assigned random state; (ii) Then a molecule is chosen 
and its position is set to be the origin by translating the 
coordinate system. In this way we build an equilibrium 
system with one molecule initially localized at the origin. 
The molecule at the origin is hereinafter referred to as 
"the tagged molecule". By assigning respectively all N 
molecules to be the tagged one in this way, i.e. setting 
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FIG. 1: The PDF p m (r, t) of the tagged molecule at time t = (a) and t = 15 (b). (For the sake of presentation, the coarse- 
grained results over a grid of squares of size 4x4 are plotted.) The intersection of the plot in (b) with plane y = is shown 
in (c) for a close look, where p m {x,y,t) = p m {r,t). (d)-(f) and (g)-(i) are the same as (a)-(c) but for the ensemble averaged 
energy density distribution p e (r,t) and the spatiotemporal correlation function c e (r, i) of the energy fluctuation respectively. 
The ensemble average for all the three cases is evaluated over 3 x 10 8 systems. 



their initial positions to be at the origin one by one, a 
subensemble of N equilibrium systems, each has a tagged 
molecule initially residing on the origin, is thus prepared, 
(iii) By repeating step (i) to have different realizations of 
the seed equilibrium system followed by (ii) to generate 
the corresponding subensemble, we then build the whole 
ensemble whose member number can be large enough (up 
to ~ 10 s ) for satisfying statistic results. It is worth not- 
ing that building the ensemble in such a way is not new; 
a similar idea was once employed by Helfand to establish 
the energy diffusion theory of fluids [4| • 

To investigate that as the system evolves, how a 
molecule (represented by the tagged molecule) diffuses 
and how the energy it carries initially spreads over the 
whole system, we will study in particular the following 
three processes with the ensemble prepared: 

A. The self-diffusion, or mass diffusion process, which 
can be accessed by calculating the PDF, denoted by 
p m (r, t), of the tagged molecule. It is defined as 



p m (r,t) = (<5[r-n(t)]), 



(2) 



where (•) denotes the ensemble average, rj(t) denotes the 
position of molecule i at time t and number one molecule 
represents the tagged molecule. As initially the tagged 
molecule is located at the origin, i.e., ri(0) = 0, we have 
p m (r,0) = <5(r). 

B. The dispersion of the energy initially carried by the 
tagged molecule. For this aim we consider the ensem- 



ble averaged energy density distribution (EAEDD) of the 
systems, i.e. 



(E(r,t)) = (H^Slr-nit)}) + <£ Hj (t)*[r ~ 

J'=2 

(3) 

Initially, as the origin is occupied exclusively by the 
tagged molecule, we have (ifi(0)<5[r — i"i(0)]) = E5(r), 
where E = (Hj) is the average energy of a molecule 
in the gas. On the other hand, as the rest area of the 
space (r ^ 0) is occupied uniformly by other N — 1 
molecules, the EAEDD they contribute to, i.e. the sec- 
ond term on the r.h.s of Eq. ([3]), equals a constant 
rj = E(N—1)/S. Hence initially (E(r, 0)) is characterized 
by a peak at the origin with a flat background. As the 
system evolves, while the portion of (E(r, t)) contributed 
by the tagged molecule may spread out from its initial 
(5-function, that by the other N — 1 molecules remains 
to be rj. Therefore we can use the reformed distribution 
p e (r,t) = ((E(r,t)) — rj)/E to capture the dispersion of 
the energy initially carried by the tagged molecule. This 
is the key technique of this work; with it the energy trans- 
ferred to others from a single molecule can thus be traced 
by checking the EAEDD of the whole system, making it 
possible to study the former conveniently. 

C. The energy fluctuation correlation. The spa- 
tiotemporal correlation function of the energy fluctua- 
tion [n| is a useful tool 11, [l2| in tracing how the 
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energy transfers 13j, [l4|. When the energy initially 
concentrated at the origin (the cause) is transferred 
to position r at time t, it will induce a solid cor- 
relation (the effect); we have to expose the correla- 
tion induced exclusively by this causality. As the gas 
model we study here is a microcanonical system with 
the total energy conserved, we have at any time r that 
Ef=i&H 3 (T) = and thus AHi(t) J2j=i &Hj( T ) = 
A#i(T)Affi(T) + J^^AHx^AH^t) = 0. Here 
AHj(t) = Hj(t) — E. Considering the ensemble aver- 
age, we then have (AHii^AHjir)) = (AiJ 2 (-r))/(A-l) 
since the gas is homogeneous, which indicates that there 
is a trivial correlation between any two molecules in- 
duced by the conservation of the energy. To get rid of it 
we consider instead (AHi(0)^2^ =1 AHj(t)S[r — ij(t)]); 
As initially (at t = 0) it has a center of (5-function 
form (Ai?i(0)A/fi(0)(5(r)) and a flat background p = 
-(A#i(0)Aiii(0))A7(7V - 1)5, we accordingly define 
the spatiotemporal correlation function as c e (r, t) = 

{AH 1 (0) Y^=i AH j W^[ r ~ r j (*)]> — M to explore the cor- 
relation at position r and time £ induced by the initial 
energy fluctuation AiZi(O). 

The main results are summarized in Fig. 1. First of 
all Fig. 1 (a)-(c) are for p m (r,t), the PDF of the tagged 
molecule. Initially it is a 8- function as seen in Fig. 1(a), 
in agreement with the fact that the tagged molecule is 
located at the origin at the beginning. Later it develops 
into a Gaussian distribution [Fig. l(b)-(c)]. As a double 
check we have also studied the squared displacement of 
the tagged molecule and found that it depends on time 
linearly after a transient time of t » 10; i.e., (|r| 2 (i)) = 
ADt with D — 7.08 ± 0.01 as suggested by the best linear 
fitting (|r| 2 (i)) over 10 < t < 80. These results are clear 
evidence that the self-diffusion of a molecule is normal in 
our system. 

Fig. 1 (d)-(f) show the energy transfer behavior. The 
(5-function seen in Fig. 1 (d) indicates that the energy 
we are interested in is initially located at the origin. 
However, in its later development, a distinctive differ- 
ence from the self-diffusion of a molecule can be identi- 
fied: Instead of a Gaussian distribution, p e (v,t) features 
a growing "crater", i.e., a ring ridge (where p e (r,t) > 0) 
moves outwards leaving behind a dip (where p e {?, t) < 0) 
in center [see Fig. 1(f)]. 

Fig. 1 (g)-(i) show the correlation function c e (r,t). It 
can be found that c e (r, t) reveals the same features of the 
energy transfer process seen in p e (r, t). The crater struc- 
ture in c e (r, t) appears at the same position and expands 
outward with the same velocity as that in p e (r,t). As 
in the center dip region of p e (r,i) we have (AHj(t)) < 
[see Fig. 1(f)], which implies (Ai7i(0)AiZj(t)) > 0, the 
center peak in c e (r,i) is consistent with the dip in p e (r,t). 
It is interesting to note that between the crater and the 
center peak there is a small region showing a negative 
correlation; This feature is different from that observed 
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FIG. 2: (a) The time dependence of the radii characterizing 
the crater structure in p e (r, t), corresponding to the top ring 
of the ridge (red bullets) and the opening of the center dip 
(blue stars) respectively. The best fittings (dashed lines) sug- 
gest their expanding speeds vh ~ l.l^s and v n ~ 0.7^. (b) 
The time dependence of the positive and negative potion of 
energy, i.e. E+ (red bullets) and _E_(blue stars), correspond- 
ing to the integrated energy density distribution over region 
r| > r v and |r| < r v , respectively. 



in the lattice models, where c e (r, t) is always positive and 
as a consequence can be employed to represent the PDF 
of the energy diffusion [l(| ■ 

Now let us have a closer look at p e (r,t), the main re- 
sult of this paper. Two key geometric parameters char- 
acterizing the crater structure of /? e (r, t) as shown in Fig. 
l(e)-(f), denoted by r v and rn respectively, are the ra- 
dius of the intersection ring on which p e (r, t) = and 
that of the ring of the ridge top where /? e (r, t) takes the 
maximum value [see Fig. 1(f)]. It is found that they both 
depend on time linearly [see Fig. 2(a)], but however, the 
speed of the ridge top, represented by v>h = dru/dt, is 
different from that of the opening of the dip v v = dr^/dt: 
The best linear fitting results suggest vh ~ 1-6^. Just 
as a comparison, it is interesting to note that the two 
speeds are comparable to the speed of sound, a macro- 
scopic characteristic; i.e., vh ~ and v n w 0.7^, 
where v s = \Jc p ksT / \jc v m is the sound speed of the 2D 
ideal gas of the same molecular mass and density and at 
the same temperature as our model. On the other hand, 
as (E(r, t)} — r\ describes how the initial energy carried by 
the tagged molecule transfers away, the interesting fact 
that p e ( r , t) has a negative center suggests that during 
this process some energy of the neighboring molecules is 
brought away in addition. This additional portion of en- 
ergy is given by — E_ where E_ = E Ji r | <r p e (r,t)dr. 
Similarly, the total positive energy carried by the bulk of 
the ridge is given by E + = E j^ >r p e (r, t)dr. Due to the 

conversation of the energy we have always E+ + E- = E. 
Fig. 2 (b) shows the time dependence of £L and that 
of E + ; initially E-(E + ) decreases (increases) but after 
a transition time it reaches a constant. In other words, 
eventually the total energy brought away by the ridge is a 
constant and is larger than the energy initially the tagged 
molecule carries. Together with the results of vh and v v , 
they suggest clearly that rather than the Gaussian diffu- 
sion, the energy carried by a molecule propagates away 
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ballistically in our gas model. 

Why the energy transfer and the molecule self-diffusion 
are so different can be actually understood within the 
framework of the random walk theory. According to this 
theory, a normal diffusion occurs if the random walker 
loses its memory of the previous states completely, oth- 
erwise an abnormal diffusion may take place instead. In 
our gas system, if we focus on the tagged molecule, we 
may find that due to its frequent collisions with oth- 
ers, its memory of the initial direction of motion suf- 
fers a quick loss. This memory loss process can be 
measured by the decay of the autocorrelation function 
A(t) = (pi(O)-pi(t)) of the tagged molecule. (Herepi(f) 
denotes the momentum of the tagged molecule.) Indeed, 
as shown in Fig. 3(a), A(t) decreases exponentially with 
time, hence the motion of the tagged molecule is equiva- 
lent to that of a random walker. 

However, the information of the initial moving di- 
rection of the tagged molecule is well remembered by 
the whole system. When the energy of the tagged 
molecule transfers to others, the memory of its ini- 
tial state may transfer to the surrounding molecules as 
well. To show this we consider the correlation function 
C(0 = SjLi(Pi(0)'Pj W) as a measure of the total mem- 
ory (note that C(t) is evaluated over the whole system). 
Dividing the momentum of a molecule, say molecule j, 
into two parts; i.e., Pj(t) = p'j(t) + p"(i) , where Pj(t) 
and Pj(t) represent respectively the momentum trans- 
ferred to molecule j from the tagged molecule and other 
molecules, we then have C(t) = 5Z^L 1 (pi(0) • p'j(t)) — 
(Pl(O)-pi(O)). This is because first (pi(O)-p^'(t)) = as 

Pj(t) is independent of Pi(0) and second X)i=iPj(*) = 
Pi(0) as the momentum pi(0) is conserved in the system. 
As a result C(t) is in fact a time- independent constant, 
suggesting that though the initial moving direction will 
be forgotten quickly by the tagged molecule itself, it will 
be well remembered by the whole system in future. 

This fact implies that the energy transfer process stud- 
ied here cannot be a Markov process. To verify that 
the memory is kept during the energy transferring pro- 
cess, we reproduce Fig. 1(e) in an alternative way: In 
preparing a gas system in our ensemble as described in 
Step (ii) , we rotate additionally the coordinate system so 
that the initial moving direction of the tagged molecule is 
along axis y. This is equivalent to considering a subset of 
Helfand's subensemble where the momentum direction of 
the tagged molecule is specified as well. Fig. 3 (b) shows 
the result, from which we can see that p e (r,t) is obvi- 
ously anisotropic and suggests clearly the initial moving 
direction of the tagged molecule. 

In summary, rather than the Brownian motion, our 
simulation study of a 2D gas at room temperature shows 
that the energy carried by a molecule would propagate 
away in a ballistic wavelike manner. Considering the en- 
semble average, the profile of the transferred energy is 



FIG. 3: (a) The autocorrelation function A(t) of the momen- 
tum of the tagged molecule decreases exponentially with time; 
(b) The same as Fig. 1(e), but the direction of the initial ve- 
locity of the tagged molecule is set to be along axis y. 



found to be characterized by a ring ridge and a dip in 
center, and both expand outwards with constant speeds 
comparable to the speed of sound. The ballistic propaga- 
tion of the energy is also confirmed by the spatiotemporal 
correlation function of the energy fluctuation. 

We emphasize that our study investigates the energy 
dispersion of a single molecule in equilibrium state, hence 
in nature the observed ballistic characteristic is distinct 
from the nonequilibrium macroscopic waves such as the 



sound wave and the heat wave 15|. For example, the 
heat wave is a macroscopic relaxing phenomenon, it may 
decay in fluids due to viscosity and heat conduction, but 
the energy density distribution in the present study does 
not decay [see Fig. 2(b)]. The heat wave can also exist 
in lattice systems [16j . but again it only survives for a 
finite time due to decaying effect. The properties of the 
energy dispersion studied here are also in clear contrast 
to those of the macroscopic heat conduction sustained by 
the temperature gradient; e.g., we have also studied the 
3D counterpart of our gas model, and both a 2D and a 3D 
hard-disc gas model, but obtained qualitatively the same 
results. However, the heat conduction may have a dra- 
matic dependence on the dimensionality in momentum- 
conserved systems 17]. In spite of these difference, the 
energy dispersion properties of a molecule must have 
underlying implications to various macroscopic energy 
transport behavior. In this respect the mode coupling 
theory of hydrodynamics, which has been shown powerful 
in bridging the microscopic and macroscopic descriptions 
of fluids 18], may provide deep insights. 

Finally we would like to suggest a possible laboratory 
testing scheme of the energy transfer process bases on the 
energy correlation function discussed (see Fig. 1 (g)-(i) 
for example simulation results). The requisite laboratory 
technique is the measurement of the simultaneous posi- 
tion and velocity of a particle immersed in a gas (fluid) 
[lit . Given this, a sample of Af/i(0) • AHj(t) can be ob- 
tained by making the measurement to a certain immersed 
particle at a time and to another after time t. Repeating 
this data collecting process till a sufficient large amount 
of samples are available; c e (r, t) can then be evaluated. 
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